#!/usr/bin/env python
# -*- coding: utf-8 -*-

from findsmoospec import findsmoospec
import pylab, os, shutil,atpy
import numpy as np
from moog import moog

from tempfile import mkstemp
from shutil import move
from os import remove, close

def replace(file, pattern, subst):
    #Create temp file
    fh, abs_path = mkstemp()
    new_file = open(abs_path,'w')
    old_file = open(file)
    for line in old_file:
        new_file.write(line.replace(pattern, subst))
    #close temp file
    new_file.close()
    close(fh)
    old_file.close()
    #Remove original file
    remove(file)
    #Move new file
    move(abs_path, file)


def spectrum_loader(specFile):
	words=specFile.split()
	ifile = open(specFile,'r')
	wave=[]; intensity=[]
	for line in ifile:
		quadruple = line.split()
		wave.append(float(quadruple[0])); intensity.append(float(quadruple[1]))
	return np.transpose(np.array([wave,intensity]))

def spectrum_loader_jacco(specFile):
	ifile = open('/home/jeffrey/PHD/ruawahia/data/jaccodata/'+specFile,'r')
	wave=[]; intensity=[]
	for line in ifile:
	    quadruple = line.split()
	    if quadruple[1]=="*********":
			continue
	    wave.append(float(quadruple[0])); intensity.append(float(quadruple[2]))
	return np.transpose(np.array([wave,intensity]))

tbl = atpy.Table('/home/jeffrey/PHD/ruawahia/data/randomtables/vL07_Marino_aBoolike.vot')

linelists="../software/MOOG2010/linelists/lin_cn42_moog_cropped.dat"
stronglines="../software/MOOG2010/linelists/strongca4226"
for i in range(6):
	pylab.clf()
	LEID=tbl.data['LEID'][i]; print LEID
	temp=int(tbl.data['temp_JS'][i]); print temp
	grav=float(int(tbl.data['logg_JS'][i]*100))/100; print grav
	FeH=tbl.data['FeH_M'][i]; print FeH
	specFile=tbl.data['specfile'][i]; print specFile
	
	j=0
	colours=["b--","r-.","g-."]
	for alpha in [0,0.10,0.29]:
		#This section creates the MOOG model. For alpha-enhanced models it modifies the model to include the right abundances
		modelatm="../software/MOOG2010/model_atm/t"+str(temp)+"g"+str(grav)+"m"+str(FeH)+"a"+str(alpha)
		os.system("atmosphy -teff "+str(temp)+" -logg "+str(grav)+" -feh "+str(FeH)+" -a "+str(alpha)+" -o "+modelatm)
		if alpha>0:
			pattern="Natoms        0          %0.1f" % FeH
			subst="""Natoms        5          %0.1f
			          8    %0.2f
			         12    %0.2f
			         14    %0.2f
			         16    %0.2f
			         20    %0.2f  """ % (FeH, 8.56+FeH+alpha, 7.58+FeH+alpha, 7.55+FeH+alpha, 7.21+FeH+alpha, 6.36+FeH+alpha)
			replace(modelatm, pattern, subst)
		spectrum=moog.run(atmosphere_model=modelatm,line_list=linelists,stronglines_list=stronglines,strong=1,wavelength_range=[4190,4260],element_number=6, XFe=0.0)
		pylab.plot(spectrum[:,0],spectrum[:,1],colours[j])
		j=j+1
	
	spectrum_smoo=spectrum_loader_jacco(specFile)
	pylab.plot(spectrum_smoo[:,0],spectrum_smoo[:,1],'k-')
	pylab.xlim([4190,4260])
	pylab.ylim([0.1,0.8])
	pylab.text(4200,0.2,str(LEID)+", "+str(temp)+", "+str(grav)+", "+str(FeH))
	pylab.legend(("[A/Fe]=0","[A/Fe]=0.10","[A/Fe]=0.29","Spectrum"),loc=4)
	pylab.savefig("figures/"+str(LEID)+"_alpha.pdf")
	#pylab.show()
	
